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In this paper we investigate the quantum phase transition from magnetic Bose glass to magnetic Bose-Einstein 
condensation induced by a magnetic field in NiCl2-4SC(NH2)2 (dichloro-tetrakis-thiourea-Nickel, or DTN), 
doped with Br (Br-DTN) or site diluted. Quantum Monte Carlo simulations for the quantum phase transition of 
the model Hamiltonian for Br-DTN, as well as for site-diluted DTN, are consistent with conventional scaling at 
the quantum critical point and with a critical exponent z verifying the prediction z — d; moreover the correlation 
length exponent is found to be v — 0.75(10), and the order parameter exponent to be ft = 0.95(10). We 
investigate the low-temperature thermodynamics at the quantum critical field of Br-DTN both numerically and 
experimentally, and extract the power-law behavior of the magnetization and of the specific heat. Our results for 
the exponents of the power laws, as well as previous results for the scaling of the critical temperature to magnetic 
ordering with the applied field, are incompatible with the conventional crossover-scaling Ansatz proposed by 
Fisher et al, [Phys. Rev. B 40, 546 (1989)], but they can all be reconciled within a phenomenological Ansatz 
in the presence of a dangerously irrelevant operator. 
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I. INTRODUCTION 

The investigation of interacting bosons in the presence of 
randomness (the so-called dirty-boson problem ) represents 
a long-standing subject in quantum condensed matter, and 
in particular a subject which remains far from being settled. 
Indeed disorder has a tremendous impact on the physics of 
bosons, introducing a significant richness in their phase dia- 
gram both at zero and at finite temperature. Intense exper- 
imental investigations on a wide variety of different physi- 
cal systems - including 4 He in porous medial cold atoms 
in disordered optical potentials 2 , and doped antiferromagnetic 
insulatorJ^, among others - have been pursued in the attempt 
to unveil the universal aspects of such a phase diagram. A 
similarly intense activity has been conducted at the theoreti- 
cal level 7 , using a variety of methods ranging from field theory 
to quantum Monte Carlo. Yet, despite more than two decades 
of intense investigation, fundamental aspects of the physics 
of dirty bosons remain obscure. Among the most salient as- 
pects under debate are the nature of the non-condensed and 
non-superfluid state induced by disorder and interactions in 
a Bose fluid - namely the compressible Bose glass or the in- 
compressible Mott glass; and the nature of the quantum phase 
transition (QPT) from a Bose glass to a strongly interacting 
superfluid Bose-Einstein condensate (BEC), induced by in- 
creasing the role of interactions while keeping the disorder 
strength fixecP. In particular, the latter transition does not 
admit a well-controlled analytical treatment, its critical expo- 
nents are generally unknown, and the existing predictions^ are 
highly debated^. 

In this paper we address the problem of the Bose- 
glass/superfluid QPT from the point of view of a model mag- 



net, NiCl 2 (i_ a; )Br2 a ;-4SC(NH2)2 (Br-doped dichloro-tetrakis- 
thiourea-Nickel, or Br-DTN), whose magnetic Hamiltonian 
can be mapped exactly onto that of 3d interacting bosons on a 
lattice with random hopping and random on-site interactions. 
A recent experiment 5 has demonstrated that this compound 
exhibits an extended magnetic Bose-glass phase in an applied 
magnetic field - such a phase is characterized by a gapless and 
magnetically disordered low-temperature state. Increasing the 
field beyond a critical value drives the system through a QPT 
from Bose glass to magnetic BEC, characterized by sponta- 
neous long-range antiferromagnetism in the plane transverse 
to the field. In particular the critical temperature to condensa- 
tion is shown to scale with the applied field as T c ~ | h — h c | ^ 
with <fi = 1.1(1), where h represents the dimensionless mag- 
netic field and h c its critical value. A seemingly puzzling fea- 
ture is that the only theoretical prediction for the <p exponent, 
stemming from Ref. 8] gives <f> > 2, in open contradiction 
with the result of Ref. 5 Using quantum Monte Carlo sim- 
ulations of the model Hamiltonian for Br-DTN, as well as of 
the Hamiltonian of DTN with a diluted magnetic lattice, we 
provide a detailed investigation of the quantum critical scal- 
ing at the Bose-glass/superfluid transition. In the case of Br- 
DTN we also investigate the low-temperature specific heat - 
both numerically and experimentally - and the magnetization 
- numerically - around the critical field. Our theoretical data 
for the zero-temperature quantum critical point (QCP) sup- 
port the picture of conventional scaling and the prediction^ of 
z = d = 3 for the dynamical critical exponent. Nonetheless 
the results for the behavior of the system at finite tempera- 
ture are seen to be globally incompatible with the predictions 
stemming from the scaling Ansatz of Ref. [8] for the free en- 
ergy, which assumes that the temperature is the only direction 
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of instability of the QCP, while disorder and interactions can 
be eliminated from the scaling. We discuss possible violations 
of conventional scaling, due to the presence of a dangerously 
irrelevant operator, and find that a fully consistent picture can 
only be obtained when assuming that such an operator affects 
the finite-temperature behavior only, preserving conventional 
scaling at T — 0. In the absence of a consistent theoretical 
treatment of the long-wavelength behavior of dirty bosons 7 , 
our results offer a first, empirical determination of the scaling 
behavior around the dirty-boson QCP. 

The structure of the paper is as follows: Se c. |H| reviews the 
main aspects of the QPT of dirty bosons; Sec.|III|shortly illus- 
trates how statically doped DTN realizes the physics of disor- 



dered lattice bosons; Sec. IV illustrates the crossover scaling 



Ansatz of Ref. [8] and its predictions for the quantum critical 
behavior; Sec. [V] discusses the numerical estimate of the dy- 
namical critical exponent and of other critical exponents at 



the dirty-boson QCP; Sec. VI presents our numerical as well 
as experimental results for the thermodynamics of Br-DTN at 



the quantum critical field; Sec. VII discusses how all of our 



results can be captured quantitatively by a generalized scaling 
Ansatz; conclusions are drawn in Sec. | Villi 



II. DIRTY-BOSON QUANTUM CRITICAL POINT 

In the absence of disorder, the physics of interacting bosons 
is generally well understood following a set of established 
paradigms. A diluted Bose gas of weakly interacting bosons 
undergoes condensation in a spatially extended state, and it is 
quantitatively described by Bogolyubov and Beliaev theory 
The QPT from the vacuum state into the condensed state, 
driven by the chemical potential, represents a Gaussian tran- 
sition in dimensions d > 1 with dynamical critical exponent 
z = 2.^ In d = 3 the onset of the critical temperature to 
condensation around the QCP follows the power-law scaling 
T c ~ |/x — A*ci (where /i c is the critical chemical potential), 
with <fi — 2/3 as predicted by mean-field theory 12 . For bosons 
on a lattice, an instability to the formation of a Mott-insulating 
state occurs for sufficiently large interactions, and the QPT 
from a superfluid condensate to a Mott insulator away from 
commensurate filling admits the same description as the QCP 
of the diluted Bose gas, while it belongs to the universality 
class of the XY model in d + 1 dimensions at commensurate 
fillingP 

The introduction of disorder to the system is accompa- 
nied by the Anderson localization of the (low-energy) single- 
particle eigenstates; such states cannot host a Bose-Einstein 
condensate in the thermodynamic limit in the presence of any 
finite interaction. As a consequence, the physics of the diluted 
Bose gas is completely altered by disorder, and a novel phase 
appears - the bosonic Anderson glass or Bose glass 8 13 - char- 
acterized by the redistribution of bosons on an extensive num- 
ber of low-lying Anderson-localized states. This fact renders 
the theoretical description of the Bose glass quite involved, 
as conventional tools for the diluted Bose gas are no longer 
available. This is a fortiori true when the interaction strength 
increases (the disorder strength being held fixed), inducing a 



QPT from the Bose glass into a strongly correlated condensed 
and superfluid state. This transition, as well as other tran- 
sitions in disordered systems^, represents an exceptionally 
hard problem in the theory of critical phenomena. Indeed the 
upper critical dimension d c , above which mean-field theory 
becomes exact, is unknown for dirty bosons; and perturba- 
tive renormalization group treatments break down 8 , given that 
they are based on the assumption of a finite d c , while d c might 
as well be infinity for dirty bosons. This means that a quantita- 
tive understanding of the behavior of interactions and disorder 
under renormalization is still lacking, as well as a rigorously 
motivated scaling form for the free energy. In the following 
we will offer evidence that the scaling Ansatz for the free en- 
ergy proposed by Ref.[8]is not appropriate. 

Another fundamental open question concerns the value of 
the dynamical critical exponent z. The derivation of the iden- 
tity z = d offered in Ref.|8]has been questioned in Ref.|9j and 
some n umeri cal calculations in d = 2 seem indeed to con- 
tradict i t ' \ although other 2d studies are instead consistent 
with it 17 """ 22 . For d = 3 much less results are available, but 
z = 3 is consistent with the existing numerical studies of the 
3d Bose-glass/superfluid transition 2223 , and, as we will see, 
also with the numerical results of this work (Sec.fV]i. 

The derivation of z = d in Ref. [8| is based on the obser- 
vation that the compressibility n = d 2 f /d^ 2 (where / is the 
free energy density and /i is the chemical potential) behaves 
as the imaginary-time stiffness of the system. An expansion 
of the singular part of the free energy upon an infinitesimal 
twist leads to the prediction that k ~ g u ( d ~ z ) t where g is the 
distance to the critical point; if the compressibility is to re- 
main finite across the transition, this would then imply that 
d = z. This derivation has been questioned by Ref. [9j on 
the basis that the main response to a twist in imaginary time 
is expected to come from the analytical rather than from the 
singular part of the free energy - this is due to the fact that, 
if the particle-hole symmetry is already broken (as it is across 
the generic Bose-glass/superfluid transition), a variation in the 
chemical potential does not alter the Hamiltonian symmetries 
in the system, and hence it should lead to negligible contribu- 
tions in the singular part of the free energy density. Therefore 
the finiteness of k across the transition should not have any 
obvious implication for the critical exponents. 

It is generally hard to believe as an accident that many nu- 
merical results appear to be consistent with z — d - although, 
as noticed in Ref. 16, only in a few cases an unbiased estimate 
of z has been provided. Without the pretension of solving the 
issue of the value of z, in the following we provide a simple 
derivation of the finite-size scaling for the compressibility at 
zero temperature, based on its relationship with the particle- 
hole gap. Indicating with E(N) the ground-state energy of 
the system having N particles, we have that the particle-hole 
gap takes the form 



A ph = E(N + 1) + E(N - 1) - 2E{N). 



(1) 



For an hypercubic system of size L d , taking the limit L,N ^> 
1 one obtains 



A 



d 2 E 



1 1 



P h 



dN 2 L d k 



(2) 



3 



where we have used the fact that fi = dE/dN, and the defi- 
nition of the compressibility k = L~ d dN/d/i. At the QCP, 
the finite-size particle-hole gap must vanish as A p h ~ L~ z , 
yielding the result k ~ L~^ d ~ z \ Hence a finite compress- 
ibility in the thermodynamic limit at the QCP would re- 
quire that z — d. The above derivation is based on the 
only assumption that conventional scaling applies (giving 
Aph ~ L~~ z ), namely that the system is below or at the up- 
per critical dimension. The scaling of k is indeed observed 
at the commensurate-incommensurate transition of the Bose- 
Hubbard model with z = 2 both in d = 1 (where the com- 
pressibility diverges 24 ) and in d — 2 (where the compressibil- 
ity jumps 27 , up to logarithmic corrections). It is also verified 
at the Mott-insulator/superfluid transition at commensurabil- 
ity (with z = 1), e.g. in d = 1 (compressibility jump at the 
Kosterlitz-Thouless transition^) and in d = 2, 3 (where the 
compressibility vanishes continuously^! 2 ^) 



m. DIRTY BOSONS AND THE PHYSICS OF DOPED DTN 

Recently a new experimental playground has become avail- 
able for the study of the physics of dirty bosons: doped mag- 
netic insulators. In such systems a Bose fluid of magnetic 
quasiparticles (corresponding to modes carrying a spin m s = 
1) can be controlled by the application of a magnetic field, 
and can be made to condense at sufficiently low temperature 
when the field exceeds a critical value 28 . The condensed phase 
corresponds to long-range (anti)ferromagnetism in the plane 
transverse to the field. Static doping of the magnetic insula- 
tor leads to disorder in the magnetic Hamiltonian, and allows 
therefore to investigate a field-induced Bose-glas s/superflu id 
transition, which has been discussed t heoreti call y 5 ! 20 ! 29 ! 30 ! as 
well as recently probed in experiments ^ 6 ! 31 ! 3 ^. 

Here we focus our attention on Br-doped DTN, which has 
been recently investigated in Ref.|5] A minimal model Hamil- 
tonian for the magnetic behavior of Br-DTN is given by 

Hbt — DTN — ^ ^ Jc,(ij) &i ' &j ~^ Jab ^ ^ &l ' &m 

(ij) c (lm) ah 

+ J2 D i(SD 2 -M B Hj2S?- (3) 

i i 

Here Sf (a = x,y,z) are S = 1 spin operators, coupled 
on a cubic lattice with antiferromagnetic interactions J c along 
the c axis and J a {, in the ab plane; D is a strong single- 
ion anisotropy. In pure DTN J c — 2.2 K, J ao = 0.18 K, 
D = 8.9 K and g = 2.26. 33 The effect of Br doping on the 
magnetic behavior can be quantitatively captured by promot- 
ing the J c couplings and the single-ion anisotropies to local, 
bond-dependent (J c ,(ij)) and site-dependent (Dj) quantities, 
related to the appearance of a Br dopant on the Cl-Cl bond 
bridging the antiferromagnetic interactions between Ni ions 
along the c axis. In particular a concentration x of Br dopants 
leads to a fraction 2x of doped Cl-Br bonds (neglecting the 
case of rare Br-Br bonds for x <C 1); Br doping is modeled as 
increasing the J c bond strength to J' c = 2.3 J c , and reducing 
the single anisotropy D to D' = D/2 on one of the two ions 



connected by the bond. In the rest of this paper we will focus 
on the doping value x = 0.075. 

A Holstein-Primakoff transformation maps the magnetic 
Hamiltonian onto a lattice boson HamiltoniarP 



rii , ,+ / ii, 

1 bib] Wl- — +h.c. 

2 3 V 2 
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y.Jijjni - l)(rij - 1) + 2J Di{rii - 1 



(H) 



-gUsH 7ij + const. 



(4) 



where bi, b\ are bosonic operators with maximum occupation 
rii = 2; here Jy = J c .(ij), Jab depending on the orientation 
of the bond. The single-ion anisotropy plays the role of an 
on-site repulsion term for the bosons (which penalizes doubly 
occupied as well as empty sites, due to the spin-inversion sym- 
metry of the Hamiltonian); the antiferromagnetic couplings 
play the role of a hopping term as well as of nearest-neighbor 
repulsion; and the role of the chemical potential is taken by 
the applied magnetic field. In particular the pure system in 
zero field is in a Mott-insulating phase with n = 1 particle per 
site (corresponding to zero magnetization on each site), and 
it exhibits a Mott-insulator/superfluid QPT at a critical field 
h c i = qhbH/ J c ~ 1.65; a second transition from superfluid 
to rt = 2 Mott insulator is reached at the upper critical field 
h C 2 = 8.69. In the presence of Br doping h cl shifts to a lower 
value h c i — 0.828(3) and the QPT connects a Bose-glass to a 
superfluid stated h c % is instead weakly affected by Br doping, 
but the corresponding transition also changes to a Bose-glass- 
to-superfluid one. 

To probe the universal character of some of our findings, 
we consider as well the case of site dilution of the magnetic 
lattice of DTN, leading to the following model Hamiltonian 

"Hsd-DTN = J c 2_j e i e j&i ' &j + Jab ^ e i^jSi • S m 

(ij)c (lm) ab 

+ dJ2^(SD 2 -9VbHJ2^S z . (5) 

i i 

The Hamiltonian parameters are now set to be homogeneous 
and equal to those of pure DTN, but the random variables t.- L 
have been introduced, taking the value with probability x 
and 1 with probability 1 — x. Here x is the lattice dilution. 
The phase diagram of the model Eq. [5] has been extensively 
studied in Ref.[35]for a dilution x = 0.15, and it exhibits two 
superfluid/Bose-glass transitions at fields h c i s» 1.94(1) and 
h C 2 ~ 8.52(2). Site dilution of the lattice of magnetic Ni 
ions in DTN can be achieved by Cd doping 36 . Similarly to the 
case of Br-DTN, a spin-boson transformation maps the spin 
Hamiltonian Eq. |5]l to that of a lattice gas of strongly inter- 
acting bosons hopping on a spatially anisotropic cubic lattice 
with site dilution. 

In the following we will mainly focus on the vicinity of 
the critical points h c i and h C 2- As sketched in Fig. [TJ, for 
d = 3 the T = QCPs are the origin of a line of critical 
temperatures T c (g) (g — \h — /i c i(2)|) f° r me transition to a 
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superfluid condensed phase. The critical temperature is follow 
the power-law dependence T c (g) ~ g^. In the pure system, 
</> = 2/3, stemming from the theory of the diluted Bose gas^, 
this result has also been established experimentally with high 
accuracy for DTN 34 . In the disordered case, both experiment 
and theory find </> = 1.1(1) for Br-DTN; this result is also con- 
firmed theoretically in site-diluted DTN. 35 It is remarkable to 
observe that other recent experiments on different compounds, 
expected to exhibit a field-induced Bose-glass/superfluid tran- 
sition, also report values of <fi which are compatible with what 
found in DTN. In particular a 4> ps 1 is compatible with the 
results on K-doped TIQ1CI3, 4 , as commented in Ref. [3T| (al- 
though Ref. 4 claims that an exponent cf> ~ 2 can be extracted 
from a fit to the experimental data, see also Ref. 132). Values 
of 4> compatible with those of doped DTN are also found in a 
recent experiment on piperazinium-Q^Cli-zBraOe- • 

This mounting theoretical and experimental evidence in fa- 
vor of <fi ps 1 for the T c line stemming from the dirty-boson 
QCP clashes with the only prediction existing in the literature, 
namely that of Ref. [8] giving <fi = vz > 2 if z = d. In the 
following section we will shortly review the scaling Ansatz on 
the free energy of dirty bosons from which the above predic- 
tion stems, as well as its further implications for measurable 
properties. 



IV. ONE-ARGUMENT CROSSOVER SCALING ANSATZ 

The onset of a finite-temperature transition at T c {h) start- 
ing from the QCPs at h&, h C 2 is a complex problem in the 
theory of critical phenomena, because the thermal transition 
and the quantum transition are completely different in nature. 
A comprehensive treatment of the two can be attempted via 
the crossover scaling theory of phase transitions^, which can 
account for the crossover between the quantum critical behav- 
ior at low and zero temperature and in the vicinity of the QCP, 
and the classical critical behavior close to the line of critical 
temperatures T c (h) (see Fig.[T|for a sketch). For a system of 
bosons in a disordered environment Ref. 8 propose the sim- 
plest possible (one-argument) crossover scaling Ansatz for the 
singular part of the free-energy density, namely 



fa(T,g)~\g 



2-ct 



F 



T 



T 



(6) 



where g = h — h&M) an d v, z and a are the critical expo- 
nent of the QCP. Here T is in units of J c /kg. Such a scal- 
ing Ansatz can be obtained via a phenomenological assump- 
tion for the behavior of the free energy under renormalization 
group (RG) 37 . Such a behavior is valid when the QCP is a 
fixed point under RG which only admits the temperature T 
and the control parameter g as directions of instability. All 
other directions in parameter space (the strength of disorder, 
the strength of interactions, etc.) are supposed to be stable 
directions, and the associated parameters can be eliminated 
from the scaling Ansatz. 

In order for the system to exhibit a finite-temperature tran- 
sition at T c , the function F(y) must exhibit a singularity at a 
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FIG. 1: (a) At a QPT, the quantum critical behavior at finite temper- 
ature manifests itself within a "fan" loosely bounded from below by 
crossover lines T co ~ g vz . The line of thermal transitions T c ~ <jr, 
and the surrounding Ginzburg region (gray-shaded area) of classical 
critical behavior, can be well decoupled from the quantum critical 
fan if <f> = vz. (b) Conjectured picture for the dirty-boson QPT: for 
dirty bosons <f> < vz, so that the hierarchy between T c (g) line and 
the putative crossover temperature T co to quantum-critical behavior 
is actually inverted. This would imply that the Ginzburg region and 
the quantum critical region touch each other. 



given value y c of its argument, so that one expects the criti- 
cal temperature to follow the scaling law T c = y c \d\ with 
scaling exponent <fi = vz. This means that, according to the 
scaling Ansatz Eq. km, the onset of the critical temperature 
is completely governed by the critical exponents v and z of 
the QPT. Assuming that z = d for the dirty-boson QCPP and 
given that v > 2/d according to the Harris criterion 3 -^, one 
obtains a scaling exponent <j) = vz > 2, as already mentioned 
previously. 

Further consequences of the scaling Ansatz Eq. [6] involve 
the thermodynamic behavior of the system along the quantum 
critical trajectory, T — > for g = 0. In fact one has that the 
specific heat C scales as 



C - T 



d 2 f s 



dT 2 



3=0 



(7) 



with 



2 - a 
x c = 1 



(8) 



while the uniform (field-induced) magnetization scales as 



m(T) — m(0) 



dg 



9 =0 



(9) 
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with 



(10) 



The scaling Ansatz Eq. |6]) is completely consistent with the 
validity of conventional scaling at the QCP, and in particular 
of hyperscaling, v{d + z) = 2 — a, so that one further obtains 
x c = (d + z)/z- 1 = landa; m = 2 - l/(uz) > 3/2. 

The prediction <f> > 2 is based on the assumption that z = d, 
which, as discussed in Sec. [II] is still a matter of debate. Hence 
one could argue that the prediction of Ref. 8 can be reconciled 
with the observations that <j> ta 1 in disordered magnets by 
assuming that z = d is indeed violated. In the next section we 
will provide an accurate study of the quantum critical scaling 
at the T = QCP for both Br-doped DTN and site-diluted 
DTN, testing the validity of the z — d prediction. 
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V. FINITE-SIZE SCALING FOR THE BR-DTN 
HAMILTONIAN 

In this section we provide a finite-size scaling analysis of 
quantum Monte Carlo results for the dirty-boson QCP exhib- 
ited by the models for doped DTN (Br-doped or site-diluted). 
Our quantum Monte Carlo results have been obtained mak- 
ing use of the Stochastic Series Expansion (SSE) approach 39 
on L 3 lattices with L = 12, 14, 16 and 18. All results are 
averaged over ~ 300 disorder realizations. The quantities of 
interest here are: 1) the correlation length along the c axis and 
in ab plane, obtained via the second-moment estimator 

_ L I S X (Q) 

where 

S^(k) = L- 3 J2 eZk ' irz ~ r3) ( s ? S 3 + s I S j) ( 12 > 

is the structure factor for the spin components transverse to the 
applied field, (...) denotes the thermal and disorder average, 
Q = (ir,Tr,ir) is the ordering vector, and q = (2n/L)c for 
p = c and q = (2ir/L)d (or (2n/L)b) forp = ab; 2) the spin 
stiffness 

».= *g<W*> (13) 

where W = (W a , Wb, W c ) is the winding number of bosonic 
worldlines along the three spatial directions^, 3) the squared 
order parameter, defined as m 2 s = S ± (Q)/L 3 . 

We have applied a /3-doubling procedure 41 , particularly fit 
for SSE, which consists in decreasing exponentially the tem- 
perature until convergence of the above quantities is reached 
towards their T = value. Fig. [2] shows that, for the system 
sizes under investigation, this approach leads successfully to 
the physical T = value of the quantities of interest for tem- 
peratures T«2* 10~ 3 J c . In particular Fig.|2]shows that con- 
vergence is achieved even in the most delicate case h = 0.83 



FIG. 2: /3-doubling plots for the correlation length, spin stiffness and 
order parameter of the Br-DTN Hamiltonian, h — 0.83. 



for Br-doped DTN, corresponding to a field very close to the 
QCP. 

If conventional scaling applies to the system - namely, if 3d 
dirty bosons are below their upper critical dimension - then 
one expects the following finite-size scaling behavior for the 
above quantities: 

i p {g,L) = LF ip [\g\L^) (14) 

Ps (g,L) = L d +*- 2 F Ps (\g\L^) (15) 

m 2 s (g,L) - L 2 ?/»F ms (VAL 1 '") (16) 

For finite-temperature data, one should also include 
the temperature dependence of the scaling functions as 
FodglL 1 /" , (3L~ Z ) for the observable O, and hence one 
would need a prior knowledge of the z exponent to be able 
to obtain the universal functions from collapse of the univer- 
sal data. On the other hand we can omit the temperature ar- 
gument of the scaling functions, having systematically elimi- 
nated thermal effects from our data; this will allow us to give 
an unbiased estimate of the critical exponents. 

A. Br-doped DTN 

We can make use of Eq. ([14} to locate the critical field 
via the intersection of the £ P /L curves plotted as a function 
of h. This is shown in Fig. [3] The crossing point for the 
iab/L curves is located at h c i = 0.830(15), with a fairly 
large error bar due to the shallow form of the curves around 
the critical field. In particular this value is in very good agree- 
ment with the value h c \ = 0.828(3) obtained by extrapola- 
tion of the T c (h) curve to zero temperature^. On the other 
hand, we observe that the crossing point for the £ C /L curves 
is at a slightly larger value, h c i — 0.850(5). We attribute 



6 



h c =0.830(I5) 



o 


L=12 


o 


L=14 


o 


L=16 


o 


L=18 



* 

♦ si 
• r 

/ o - 



JJLP 



h c =0.850(5) 



• oo 



"■!, 



H 

i v = 0.72(7) 



0.78 0.8 0.82 0.84 0.86 0.88 




0.82 0.84 86 0.88 0.9 
h 



1 1 1 1 1 1 
- v = 0.66(5) 


i 1 i 




A 






"1,1,1,1,1, 





\h-h ,1 L 

cl 



IN 



mate. We conclude this analysis by examining the behavior 



!>J 200 
Cl 



--3.0(1) 



TT 7 " 



,i A 



,t*o° ° 0< 





12 


o o L= 


14 


o o L= 


16 


»--*L= 


18 




0.78 0.8 0.82 0.84 0.86 0.88 0.9 



2 
IN 



\h-h ,1 L 



FIG. 5: Scaling plots of the QMC data for the spin stiffness of the 
Br-DTN model. 



FIG. 3: Scaling plots of the QMC data for the correlation length of 
the Br-DTN model. 



this discrepancy to the spatial anisotropy present in the Br- 
DTN model in terms of disorder effects, given that bond ran- 
domness only affects the J c couplings and not the J a b ones. 
This asymmetry might hinder the development of correlations 
along the c axis, and hence alter the finite-size estimates of the 
critical field. Nonetheless a very good collapse is obtained in 
a scaling plot as a function of gL x l v with v = 0.72(7) for the 
£ Q & curves and v = 0.66(5) for the £ c curves (only for positive 
g for the latter, again showing an asymmetry between the two 
spatial directions). 



of p s . Imposing that the p s L d+z ~ 2 curves cross at the esti- 
mated h c i, we obtain that z = 3.0(1), fully consistent with 
z = d (see Fig. [5}. A full scaling plot of p s with v ps 0.7 does 
not give a very satisfactory collapse far from the critical point. 
We attribute this lack of scaling to the asymmetry between 
the ab plane and the c axis: indeed p s contains the winding 
number fluctuations along all spatial directions, whose scaling 
properties for the sizes we have examined might be different. 
We postpone to future work a detailed analysis of the winding 
number fluctuations along specific spatial directions. Another 
possibility is that corrections to scaling must be taken into ac- 
count in the behavior of p s for the system sizes considered 
here (although they did not seem necessary to account for the 
scaling of the other observables). 
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FIG. 4: Scaling plots of the QMC data for the squared order param- 
eter of the Br-DTN model. 



FIG. 6: Scaling plot of the QMC data for the spin stiffness of the 
Br-DTN model using z = 1.6. The scaled curves do not appear to 
cross at the same field. 



To better refine our estimate of h c \ we consider the scal- 
ing of the order parameter, depicted in Fig. |4] We observe 
an excellent collapse for h c \ = 0.830(5), v = 0.70(5) and 
f3 = 0.95(5). On the other hand, we cannot obtain such a 
good collapse both for negative and positive g values by us- 
ing the larger /i cl estimate given by the scaling of £ c . Hence 
in the following we will make use of the h c i ~ 0.83 esti- 



Relying on the conclusions of Sec. |IV| one would argue 
that vz = 4> would rather imply that z = <f>/v « 1.6 (us- 
ing (f> = 1.1 and v = 0.7). On the other hand, when using 
this value for z we obtain that the p s (L)L d+z ~ 2 curves do not 
cross at one and the same critical field - see Fig. [6] Willing 
to identify an estimate of the critical field from the crossing 
of the two largest system sizes (L — 16 and 18) we obtain 
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h c i = 0.87(1), clearly inconsistent with the previous esti- 
mates. Hence the equality <f> — vz is definitely incompatible 
with our results. 



B. Site-diluted DTN 

We can repeat a similar analysis to that of Br-doped DTN 
for QMC results on site-diluted DTN with x = 0.15 dilution. 
All results are summarized in the scaling plots of Fig. [7] re- 
ferring to the case of both h c \ and h C 2- The central result 
is that all critical exponents appear to be consistent between 
the two transitions. In the case of site dilution, the asymme- 
try between the ab plane and c axis in terms of the disorder 
distribution is absent, and indeed we do not observe a signif- 
icant discrepancy in the critical fields estimated via and 
£ c ; moreover we obtain a good collapse of the p s data, unlike 
what observed in Br-DTN. 

Even more importantly, the critical exponents estimated for 
site-diluted DTN are fully consistent with the estimates ob- 
tained above for Br-doped DTN, as well as with the esti- 
mates of Ref.|23]for the Bose-Hubbard model with a random 
chemical potential. The latter reference gives v = 0.70(12), 
z = 3, and r\ w —1, from which we can extract j3 = 
v(d + z — 2 + rf)/2 ~ 1.05. Hence we can conclude that 
none of the results obtained so far are specific of the kind of 
disorder we have considered, and that they rather reflect the 
universal behavior of interacting bosons in d = 3 in the pres- 
ence of short-range correlated disorder. 



VI. THERMODYNAMICS ALONG THE QUANTUM 
CRITICAL TRAJECTORY OF BR-DTN: NUMERICS AND 
EXPERIMENT 

A further test of the crossover scaling Ansatz comes from 
the direct determination of the exponents xc and x m from 
numerical simulations and experiments. In the following we 
describe the extraction of these exponents from our quantum 
Monte Carlo simulations close to the h c \ critical field for the 
Br-DTN Hamiltonian, as well as from specific heat experi- 
ments on Br-DTN. 



A. Quantum Monte Carlo 

In Fig. [8] we analyze the scaling of the energy per spin, 
E(T), as obtained from our QMC simulations. The ther- 
mal energy is expected to scale as a power law close to the 
QCP, E(T) = £7(0) + A E T XE where x E = x c + 1. We 
fit our QMC data for two field values close to the estimated 
hci = 0.828(3) T value, h cl - 0.04 T and h cl + 0.007 T, 
and for two system sizes, L = 14 and L = 18. We make 
use a windowing technique, extracting the fitting parameters 
over a decreasing temperature window [0, T] in order to ob- 
tain the asymptotic quantum critical behavior for T — > 0. The 
results of this fitting procedure are shown in Fig. [8] For both 
system sizes and field values we obtain that the x E estimate 



converges to a value xe — 2.8(2) as the temperature window 
for the fit is progressively reduced, whence we extract the es- 
timate x c = 1.8(2). 

A similar procedure is used to analyze the magnetization 
m(T), which is fitted to the form m(T) = m(0) + A m T Xm 
where A m is a constant. The result of the windowing analysis 
is shown in Fig.|9j yielding an exponent x m — 1.8(1). 



B. Specific heat measurements 

We can perform a similar analysis to the experimental spe- 
cific heat data of Br-doped DTN close to the experimental H c i 
value, H cl = 1.07(1), namely for H = 1, 1.1 and 1.2 T. The 
specific heat has been measured in a dilution refrigerator us- 
ing the quasi-adiabatic heat pulse method. The analysis of the 
low-T experimental data for the specific heat is complicated 
by the appearance of a nuclear Schottky anomaly, which intro- 
duces a term that scales roughly as H 2 /T 2 in the specific heat 
- at least at sufficiently small fields - as already observed in 
pure DTNp^ElThe Schottky anomaly contribution to the spe- 
cific heat at H c2 swamps the low-T electronic contribution, 
and thus we focus on data near H c i . 

We fit the specific heat data to the form 



C(T) 



a n H 
T 



[l - tanh 2 (b n H/T)] + a e T xc (17) 



where the first term mimics the temperature behavior of the 
specific heat of paramagnetic nuclear spins, and the second 
term is the electron-spin contribution. In particular, we fit the 
specific heat data at H = 1.1 T over a variable temperature 
range [0,T] using the same windowing technique as for the 
theoretical data (see Fig.[T0fb)). This gives us stable fit param- 
eters as T is lowered, and for T = 0.6 K we obtain a very good 
fit with a n = 0.0155(2) [JK/^mol)] 1 / 2 , b n = 0.034(2) 
K/T, a e = 1.004(7) J/(K 1+XC mol) and x c = 2.00(1) (see 
Fig.pta)). 



The other data sets are instead fit differently: we use as fit- 
ting parameters a e and xc only, while a n and b n are kept fixed 
- as the nuclear contribution depends explicitly on the field. 
For if = 1.0 T we obtain a e = 1.03(1) J/(K 1+XC mol) and 
x c = 2.15(1), while for H = 1.2 T we find a e = 1.004(1) 
J/(K 1+XC mol) and x c = 1.89(1) - see Fig. [Tofc-d). For 
both fields the nuclear part is quantitatively reproduced by us- 
ing the explicit field dependence of Eq. ( fTTj ). In particular for 
H = 1.2 T the experimental data deviate significantly from 
the fitting form for T m 0.2 K, due to the apparition of an 
anomaly associated with the onset of long-range order in the 
electronic spins (the value of T c sa 0.2 K at H = 1.2 T is 
in good agreement with the T C (H) curve reported in Ref. O. 
Nonetheless above T c the specific heat appears to follow a 
power-law scaling characteristic of the quantum-critical be- 
havior. 
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FIG. 7: Scaling plots for the site-diluted DTN Hamiltonian with dilution x == 0.15. 



C. Discussion 



VII. GENERALIZED SCALING ANSATZES 



The specific heat exhibits a power-law behavior along the 
quantum critical trajectory with an exponent xc — 1.8 — 2 
determined both by numerical simulations and by experiments 
on Br-DTN. Such an exponent is a clear signature of the quan- 
tum critical behavior, as it clearly differs from the power-law 
behavior with xc = 3 appearing in the BEC phase 5 . More- 
over, this result is inconsistent with the prediction xc = 1 
coming from the scaling Ansatz Eqs. |6| and Eq. ([8]) with 
z = d. The magnetization also exhibits a distinct quantum 
critical behavior with exponent x m a 1.8, Such a behavior is 
markedly different from that exhibited in the BEC regime, in 
which the thermal magnetization initially decreases as — T 3 / 2 
when T grows above zero^l In turn the quantum critical expo- 
nent x m found with QMC is not consistent with the prediction 
x m fa 1.5 from Eqs. (|6| and Eq. ( fT0| ) with z = d and v pa 0.7. 
These results, together with the value <f> = 1.1(1), call for a 
revision of the scaling behavior of the free energy governing 
the region around the QCP. 



As already mentioned in Sec.|IIJ the long-wavelength effec- 
tive action describing the physics of dirty bosons around the 
QCP does not lend itself to a well controlled RG analysis, due 
to the non-perturbative role of disorder, and to the fact that the 
upper critical dimension appears to be d c = oo.™ Therefore 
an educated guess of the scaling Ansatz for the free energy is 
not straightforward. In the following we discuss two attempts 
at a generalization of the scaling Ansatz, Eq. based on 
the possible presence of a dangerously irrelevant term in the 
effective action, namely a direction in parameter space which 
cannot be simply eliminated in the RG flow. We cannot say a 
priori whether the term in question is represented by the inter- 
action or by the disorder, or by a combination of the two. The 
generalization of the scaling Ansatz, although introduced in 
a heuristic way, can lead to relationships between the critical 
exponents which are compatible with the numerical and ex- 
perimental values cited above for Br-DTN. A fully consistent 
picture with our results can only be obtained if a dangerously 
irrelevant operator enters in a weak way in the scaling, only 
affecting the finite-temperature behavior. 
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FIG. 8: (a)-(b) Thermal energy E(T) - £(0) of the Br-DTN Hamil- 
tonian close to the lower critical field h c \. The dashed lines cor- 
respond to a fit with xe = 2.85; (c)-(d) xe exponent obtained by 
fitting the data in (a)-(b) over a temperature window [0, T]. 
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FIG. 9: Upper panel: thermal magnetization of the Br-DTN Hamilto- 
nian close to the lower critical field h c \ ; the dashed line corresponds 
to a fit with x m — 1.8. Lower panel: x m exponent obtained by 
fitting the m(T) data over a temperature window [0, T]. 



A. Scaling Ansatz in presence of a dangerously irrelevant term 

We consider the free energy density f s (g, T, u), depending 
on the distance g to the critical point, the temperature T and 
the coefficient of a dangerously irrelevant term u. Close to 
the QCP, f s is assumed to scale under m cycles of an RG 
transformation - scaling the spatial dimensions by a factor b 



FIG. 10: Specific heat of Br-DTN close to the lower critical field 
H c \. The solid lines are fits as described in the text. 



and the imaginary-time dimension by a factor b z - as 

f s (g, T, u) ps b- m{d+z) f s (g b my <> ,T b rnyT ,u b~ my ") . 

(18) 

where y g , yx and y u > 0. Stopping the RG flow when 
gb mVs > ss 1, namely when the correlation length is renormal- 
ized to order 1, £/b m ~ 1 with y g = 1/v, and identifying 
yx = z, one obtains the standard scaling Ansatz 



f s (g,T,u)^g^ d+z ^^,ug^ 



(19) 



Letting u — > leads to singularities in f s or in its derivatives. 
Following Ref.04]one can assume that $(2/1,2/2) has a singu- 
larity when 2/2 — > 0, in the form 



$(yi,?/2) w A{y iy P) y 2 ■ 



(20) 



where p and /i are unknown exponents, except for the fact 
that /i > 0. This assumption immediately leads to the scaling 
Ansatz 



f s (g,T,u) w g"~ a A 



T 

If 



(21) 



where 2 — a = v(d + z) — vy u [i < v(d + z) (in possible 
violation of hyperscaling) and <fr = vz — py u v < vz (if p > 
0). Analogously to what discussed in Sec. |IV| the thermal 
critical line must correspond to a singularity in A(y) at y c , 
so that T c ~ g® . The inequality <\> < vz is clearly verified 
by our numerical and experimental values for the cf>, v and z 
exponents of Br-DTN. 
Rewriting Eq. pi) as 



2-c« 9(7.- a) 

f s (g,T,u)^T^u^B 



/4> u p/<t> 



(22) 



and applying Eqs. |7| and (pj, one obtains the following ex- 
ponents for the along the quantum critical trajectory 



1 — a 



x c = av. 



(23) 
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Due to the violation of hyperscaling, a is a priori unknown. 
One can determine a by imposing that x m — 1.8(1), obtain- 
ing a — —0.98(1), and hence xq = 1-7(3), in good agree- 
ment with the value obtained numerically and experimentally 
forBr-DTN. 

Yet the above result has the major drawback of requiring a 
strong violation of hyperscaling at the QPT, 2 — a < v(z + d). 
Violation of hyperscaling is in open contradiction with the 
analysis of Sec. |VJ which is successfully based on conven- 
tional scaling forms for a QPT below the upper critical dimen- 
sion. This suggests that a possible dangerously irrelevant op- 
erator enters in a more subtle way, leaving the quantum critical 
behavior unaffected and only manifesting itself in the finite- 
temperature singularity of the free energy. 



B. Continentino's scaling 

The previous scaling Ansatz, albeit giving the prediction 
Eq. ( f23| l which is fully compatible the exponents for Br-DTN, 
is not the most general Ansatz in presence of a dangerously 
irrelevant term (as also pointed out in Ref. l44t . In the ab- 
sence of disorder, the boson-boson interaction is known rig- 
orously to play the role of a dangerously irrelevant terrrPH and 
yet Eq. ( p3| is not satisfied: indeed one has that = 2/3 (see 
Ref.[l2l and a — (Gaussian fixed point), leading to xc = 2, 
while the exact result is xg = 3/2, as also measured experi- 
mentally in undoped DTN! 4 -^ 

ContinentincP^ has proposed an alternative scaling Ansatz 
to describe the singular part of the free energy close to the 
classical critical line (namely in the so-called Ginzburg re- 
gion, see Fig. [T| in the presence of a dangerously irrelevant 
term: 

f s{g ^u)^m 2 - a ^{j^y-)j (24) 

where g(T) — g — uT 1 ^ is the distance to the thermal crit- 
ical value of the control parameter, shifted with respect to 
its T = value g. The function ^(t) is supposed to have 
the following properties: 1) ^>(t — > 0) = const., so that the 
Ansatz reproduces the conventional quantum critical behav- 
ior f s ~ g 2 ~ a at T = 0; 2) *(f ->■ oo) ~ t x , where 
x = (cut — a)/(vz), and ar is the critical exponent of the 
thermal transition, so that f s ~ \g(T)\ 2 ~ aT when approach- 
ing the thermal critical line. In particular in the above Ansatz 
the operator u disappears from the scaling when T — > 0, leav- 
ing the QPT unaffected. Hence hyperscaling can be safely 
assumed at T = 0. 

In general Eq. p4| is only guaranteed to work in the 
Ginzburg region 46 , but one can test the consequences of its 
extension to a broader range of validity, including the quan- 
tum critical trajectory, as done in Ref. [45] Using hyperscaling 
at both the classical and quantum phase transition, this leads 
to the predictions 

x c = v T d ( \ - — ) + - (25) 

\<p vz J z 



and 

x m = x c + 1 - \ (26) 

where v? is the critical exponent at the thermal transition. 

Using (j) = 1.1(1), v = 0.70(5], z = d and v T = 0.669(1) 
for the 3DXY universality class 47 , we obtain xc = 1.9(2) 
and x m — 2.0(2), fully compatible with our numerically and 
experimentally determined values. 

A possible argument can be formulated to justify the ap- 
parent success of the above Ansatz a posteriori. This Ansatz 
implies the rather unconventional situation, Eq. ( |2"5| l, in which 
both the quantum and the classical critical exponents enter in 
the scaling along the quantum critical trajectory 46 . This is 
not what one would typically expect, given that the classi- 
cal critical exponents pertain to the Ginzburg region, while 
the quantum critical exponent pertain to the quantum critical 
fan above the QCP The Ginzburg region surrounds the clas- 
sical critical line T c ~ \g\ , while the quantum critical fan 
is (loosely) lower bounded by a crossover line T co ~ g vz re- 
flecting the typical energy scale of the spectral features which 
vanish at the QCP. If cj) = vz, as predicted by conventional 
scaling, Eq. d5J, the two above regions can be well separated 
(Fig. [TJ a ))- On the other hand, if tf> < vz (as in the dirty- 
boson case, in which vz s» 2(f)) the quantum critical regime 
is not lower-bounded by T co but rather by T c , and hence 
the quantum-critical region and the Ginzburg region are com- 
pletely contiguous, so that one could imagine that the quantum 
critical trajectory could be influenced by the classical criti- 
cal behavior (Fig. |TJb)). This very unconventional scenario 
is purely speculative at this stage, and it will require more 
investigations^. 

C. Alternative scaling 

To avoid the shortcomings of Continentino's scaling - 
namely the extension to the quantum-critical trajectory of a 
scaling form expected to be valid close to the classical critical 
line only - we can propose a further scaling Ansatz which gen- 
eralizes Continentino's one while preserving its virtues. In- 
deed Eq.|24]is consistent with the more general, two-argument 
scaling Ansatz 

/,(T, 5 ,u)«^*>0( 5 ^,uT*-£) . (27) 

To obtain hyperscaling at T = 0, we simply need to assume 
that G{yi,yt) — > const, when y\ — > oo and yi — > 0. 
We rewrite the above scaling form as 

f s (T,g,u)^T^ ^_0_ uT *-£) (28) 

where T -> (g/T^^y^+^Q when approaching the QCP 
at very low temperatures, namely T —> and g — > but 

g / T x/{uz) qq (implying t < T co - g uz ). 

On the other hand, when approaching the quantum-critical 
trajectory, T -t 0, g -> 0, but g/T 1 /^ -t (namely 
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T ^> T co ~ 9 UZ ), we assume a rather different form for the T 
function, namely 

Hvu 2/2) ~ vlihyi + hy-2 + ....) (29) 

when 2/1,2/2 — (/1 and /2 are constants j^l. Upon this as- 
sumption we readily obtain that, along the quantum-critical 
trajectory y\ = 0, 2/2 —> 

<f+z , x /l 1 \ 

= + + 1 t 1 (30) 

z \<f> vz J 

and 

£ m = £c - T + 1 • (31) 

9 

The 6* exponent can be determined by Eq. |30|, 9 w 0.9, using 
the numerical values of all the other exponents. Hence we can 
obtain a picture which is globally consistent with all of our 
results. 



VIII. CONCLUSIONS 

We have investigated, both theoretically and experimen- 
tally, the field-induced quantum phase transition in the mag- 
netic Hamiltonian of DTN in the presence of bond/anisotropy 
disorder or of site dilution. This quantum phase transition rep- 
resents a well controlled realization of the dirty-boson quan- 
tum critical point (QCP) in three dimensions. We have pro- 
vided compelling evidence that conventional scaling is obeyed 
at the dirty-boson QCP, but that unconventional scaling (in- 
cluding the possible presence of a dangerously irrelevant op- 
erator) is necessary to account for the thermodynamics along 



the quantum critical trajectory, as well as for the power-law 
dependence of the line of critical temperatures on the applied 
field. Our results are compatible with the prediction z = d, 
but incompatible with the prediction <f> = vz for the expo- 
nent of the power-law scaling of the cr itical temperature, as 
also found in recent experiments ^ 6 ! 31 ! 3 ^. Within a generalized 
scaling Ansatz plays the role of an independent exponent, 
whose value, combined with that of quantum critical expo- 
nents as well as classical critical ones, can be quantitatively 
related to the exponents governing the power-law scaling of 
the specific heat and of the magnetization along the quan- 
tum critical trajectory. A recent experiment on piperazinium- 
Cu2(Cli_ K Br a; )6 has measured the order-parameter exponent 
via neutron scattering at low temperature, 6 giving (3 rj 0.5. 
Our evidence is that this exponent is not compatible with the 
dirty-boson universality class in d — 3 (see also Ref. l23l : 
more investigations will be necessary in this direction. 

While a conclusive picture on the dirty-boson quantum crit- 
ical point is still far from obvious, our results provide a quan- 
titative scenario which a well-controlled theory of the dirty- 
boson effective action 8 should be able to reproduce. More 
generally, our results further confirm that disordered quantum 
magnets represent an invaluable tool for the investigation of 
complex phenomena of interacting bosons, and that they stand 
among the best candidate systems to unveil the universal fea- 
tures of dirty bosons. 
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